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M^TODO DE VISUALIZACI6N DE FLUJOS MULTIFASICOS USANDO 
TOMOGRAFfA DE CAPACITANCIA ELECTRICA 

DESCRIPCldN 

CAMPO T^CNICO DE LA INVENCI6N 

La presente invencion se relaciona con un nuevo metodo de reconstruccion 
de imagenes para la visualizaci6n de flujos multifasicos usando tomografia de 
capacitancia electrica (TCE), para emplearse en medidores tomogr^ficos de flujo 
multifasico en la cuantlficaclon de los diferentes fiuidos (como gas, acelte y agua) 
producidos por un pozo petrolffero. Tambien se usa para la visuallzaci6n y 
monitoreo de otros flujos y procesos multifasicos que ocun'en en la industria tales 
como: optimizaci6n del diseno y operaci6n de separadores, optimizacl6n del 
disefio de unldades de desintegracidn de fiuidizaci6n catalrtica (FCC) y sistemas 
de lecho fluidizado, y optimizacion del diseno de reactores de lecho fijo. 

ANTECEDENTES DE LA iNVENCl6N 

En t6rminos generates, la tomografia sirve para obtener una Imagen del 
corte transversal de un objeto en un piano detemilnado. La tomografia de rayos X 
fue la primera en ser desarrollada (en la decada de 1970) y actualmente se usa 
de manera rutinaria en el area medica, asi como en algunas aplicaclones 
industriales (por ejempio, inspeccion Interna de piezas y deteccion de fallas en 
materlales). 

Posterlormente, se han venido desarrollando nuevos metodos de tomografia 
orientados m^s Wen a procesos Industriales, conocidos con el nombre generico 
de tomografia de procesos (Williams y Beck, 1995). La finalidad de estos 
metodos, que comenzaron a evolucionar a mediados de la decada de 1980, es 
obtener una imagen de la distrlbucion espaclal de las fases o componentes en un 
proceso industrial, usando Onicamente sensores extemos y sin causar ninguna 
perturbacion en el proceso, tal como se indica en la figura 1, que muestra un 
sistema de tomograffa de procesos compuesto de un tanque o tuberfa (A), un 
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sistema de adquisicion de datos (B) y una computadora de reconstruccion de 
im^genes (C). En otras palabras, la tomografra de procesos proporciona una 
manera de 'mirar' desde afuera al interior de los mismos, sin necesidad de 
alterartos ffsicamente, y esto representa un modo radicalmente distinto de obtener 
informacion estmctural del proceso a nivel global, a diferencia de los metodos 
tradicionales, que se basan en el muestreo local en un detemiinado nOmero de 
puntos. El proceso puede ocurrir en un reactor o mezcladora, en un lecho 
fluidizado, en el interior de un separador, o dentro de una tuberfa transportando 
flujo multifasico, por poner algunos ejemplos. 

Existe toda una gama de principios y tecnicas que pueden s& explotados en 
la tomografia de procesos, incluyendo metodos electricos, basados en la medicion 
de impedancia, ultrasdnlcos, de resonancia magnetica. opticos, y a base de 
radiaclon ionizante (rayos X o gama). En general, los m6todos de radlacion 
ionizante producen las im^genes de mejor definlci6n, pero son relativamente 
lentos- Por otro lado, los metodos electricos producen Imagenes de baja 
resoluci6n pero son mucho m^s r^pidos, de construcci6n robusta y de costo 
relativamente bajo. De ahf la necesidad de contar con mejores metodos de 
reconstruccion de imagenes que, tales como los de la presente Invenclon, 
permitan obtener imagenes con una mayor definicion y fidelidad. 

Particularmente en el caso de la tomografia de impedancia electrica o, 
simplemente, tomografia electrica, se ha dado un progreso muy notable en los 
ultimos anos. Este tipo de tomografia tiene dos modalidades principales: la 
tomografia de capacitancia y la tomografia de resistencia. En un sistema de 
tomografia de capacitancia (Beck et aL, 1997; Gamio, 1997; Plaskowski et aL, 
1995), (tal como se Indica en la figura 2, fonmado por un sensor (A), un sistema de 
adquisicion de datos (B) y una computadora de reconstruccion de imagenes (C)), 
nomialmente usado con mezclas donde la fase continua es aislante, se utiliza un 
sensor, tal como se indican en las figures 3 y 4, fomnado por un arreglo circular de 
electrodes (B) distribuldos alrededor de la secclon transversal a examinar, y se 
mide la capacitancia el6ctrica de todos los posibles pares de electrodos. Por 
medio de una computadora y un algoritmo de reconstruccion de imagenes 
apropiado, los datos medldos son usados para producir un mapa (o sea, una 
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imagen) de las variaciones de la constante dielSctrica (o permitlvidad relativa) 
dentro del sensor, proporcionando asf una indicacl6n grafica de la distribucion 
interna de los diferentes componentes de la mezcla. Con relacion a las figuras 3 y 
4, los electrodes (B) pueden estar colocados en el exterior de un tubo aislante (A), 

5 a fin de simplificar la fabricacidn y evitar el contacto directo con los fluidos (E) del 
proceso, Un segundo tubo metalico externo electricamente aterrizado (C) sirve 
como blindaje o pantaila el6ctrica y para proporcionar resistencia mecanlca. El 
sensor tambien cuenta con dos electrodes de guarda cilmdricos (D) que se 
mantlenen asimismo aterrizados. 

10 La patente britanica GB 2 214 640, del 6 de septiembre de 1989, describe un 

sistema de tomografia de capacitancia electrica en ei cual el metodo empleado 
para reallzar la reconstruccl6n de ImSgenes es un algoritmo de proyeccion inversa 
lineal, conocido en Ingles como 'linear back-projection' (o 'LBP'). Sin embargo, 
dicho metodo de reconstruccion produce imdigenes de calidad reiatlvamente baja. 

15 Los metodos de reconstruccion de imSgenes propuestos en la presente invenci6n 
permiten superar esta seria limitacion. 

La tomograffa de resistencia, por otro lado, esta orientada a mezclas donde 
la fase contlnua es conductora (Piaskowski et aL, 1995; Williams y Beck, 1995). 
En este case, los electrodes se colocan al ras de la superficie interior del 

20 recipiente que contiene la mezcla, y en contacto directo con los fluidos de la 
misma. Se aplican una serie de corrientes de excitacion y se miden los voltajes 
producidos, y a partir de ellos se encuentra la distribucion de conductlvidad dentro 
del sensor, la cual indica la distribucion de los diferentes componentes. 

En principio, la tomograffa de capacitancia electrica (TCE) tiene importantes 

25 aplicaciones en la medicion de flujo multifasico, y en particular de flujo bifasico 
gas-aceite, tal como se presenta comunmente a la saiida de muchos pozos 
petrolfferos. La manera tradicional de cuantificar los diferentes fluidos producidos 
por un pozo consiste en separar la mezcla de los mismos por gravedad usando 
grandes tanques, para posterlormente medir cada componente en forma 

30 individual usando medidores de flujo homogeneo convencionales. En las Oltimas 
dos decadas han surgido medidores de flujo multifasico que permiten cuantificar 
la produccidn de un pozo sin necesidad de separar la mezcla (Thorn etal., 1997). 
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Sin embargo, los medidores multifasicos disponibles actualmente sufren de una 
sensibilidad a cambios en el patr6n de flujo, a menos que se utilicen con 
dispositivos mezcladores o acondicionadores de flujo, lo cual introduce pSrdidas 
irrecuperables de presion (que finalmente se traducen en pSrdidas de energfa). La 

5 limitacion anterior se podria superar mediante la TCE, ya que esta permlte 
conocer el patr6n de flujo y ese conocimiento se podrfa utilizar para compensar la 
respuesta de los medidores multifasicos convencionales, o, alternativamente, 
tambien es posible en principio disenar un nuevo tipo de medidor multifasico 
tomografico, basado en el analisis de las imagenes de TCE obtenidas en dos 

10 pianos transversales de la tuberia proximos entre si (Hammer et al., 1997; 
Plaskowski et aL, 1995). Adicionalmente, la TCE tiene apllcaciones potenciales en 
la visualizacion, monitoreo y control de numerosos procesos multifasicos 
industriales. 

Sin embargo, el factor principal que hasta ahora ha limitado la apllcaclon 

15 practice de la TCE es la poca fidelidad o exactitud de las imagenes obtenidas 
usando los mStodos de reconstrucci6n de imagenes disponibles, de ahi la 
Importancia de contar con m^todos mejorados como el de la presente Invencion 
(Yang y Peng, 2003). Los metodos directos simples como el de proyeccion 
Inverse lineal (LBP por sus siglas en ingles) producen Imagenes relativamente 

20 deficientes que solo proporcionan una indicacion cualitativa de la distribucion de 
componentes dentro del sensor. Tal es el case de los metodos que se describen 
en la patente britanica GB 2214640 del 6 de septiembre de 1989. 

Por otro lado, los metodos mas sofistlcados, basados en tecnicas iterativas 
de optimizacion local, por lo general necesitan de uno o varios parametros de 

25 regularizacion cuyo valor optimo depende precisamente de la imagen que se 
pretende reconstruir (y que no se conoce), ademas de que la regularizacion 
Introducida tiene el efecto de suavizar los contomos de la imagen, haciendola mas 
difusa. Tal es el case de los metodos que se describen en la patente britanica GB 
2 329 476 del 24 de marzo de 1999. 

30 Existe pues la necesidad urgente de contar con metodos de reconstruccion 

mejores y con un menor error. Como un ejempio de ello, recientemente fue 
patentado un metodo de reconstrucci6n de imagenes basado en la apllcaci6n de 
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redes neuronales artlficiales (patente norteamericana US 6,577,700 B1 del 10 de 
junio de 2003). 

En la presente invencidn, se desciiben nuevos metodos de reconstnjcci6n de 
imSgenes basados en recocido simulado y algoritmos geneticos. De acuerdo con 
las leyes de la f(sica, en especial de la electrostdtica, el sensor usado en TCE (tal 
como se observa en la figura 4) puede ser considerado como un caso particular 
de un sistema de conductores en un medio dielectrico, cuya teoria fue 
desarrollada originaimente por J. C. Maxwell (1873). En nuestro caso, los 
electrodes detectores (B) actuan como los conductores mlentras que el tubo 
aislante (A) y el contenido del sensor actuan como el medio dielectrico. Para un 
sensor de n electrodos, la carga electrica qi Inducida en los electrodes esta 
relacionadas con el potencial el^ctrico v; de los mismos por medio del siguiente 
conjunto de ecuaclones lineales simult^neas: 



16 



(1) 



9n = 



C„lVi + C^2^2 + 



C V 



20 



donde los cu son los coeficientes de capacitancia propia (o simplemente la 
capacitancia propia) del electrode i, mlentras los demas, c,y , con / ^ j\ son los 
coeficientes de capacitancia mutua (o simplemente las capacitancias mutuas) 
entre los electrodes i y j\ Puesta en forma de matrices, la ecuaci6n (1) queda 
como: 







^2 









^11 ^12 
^21 ^22 



-In 



(2) 



25 



q = Cv 



(3) 



Dado que las capacitancias mutuas tienen la propiedad de reciprocidad, esto 
as, cy = Cji, solo hay m = in(n-l) capacitancias mutuas independtentes, 
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correspondientes a la matriz triangular inferior o superior de C, y correspondlentes 
tambien a cada uno de los m pares distintos de electrodos que es posible formar 
en el sensor. 

El valor de las capacitancias mutuas es una funcidn no lineal muy 
complicada de la geometria del sistema de conductores y de la distribucion 
espacial de la cxDnstante dielectrica o permitividad relatlva e (que llamaremos 
tambien simplemente 'permitividad') de los materiales dielectricos que los 
separan. En el case del sensor de TCE, la geometria de los electrodos, la del tubo 
aislante y el valor de su constante dielectrica, son todos fijos. Por lo tanto, 
podemos decir que las capacitancias mutuas son una funcion no lineal s61o de la 
dlstribucl6n espacial de la constante dielectrica en el interior del sensor, e(x;y). Al 
problema de calcular las capacitancias mutuas correspondientes a una 
distribucion especifica de permitividad dentro del sensor se le llama problema 
directo. 

El empleo de los electrodos de guarda cilindricos (D) mostrados en la figura 
3, junto con la suposicion de que la distribucion espacial de la permitividad cambia 
poco en el sentido axial, permite representar el sensor por medio de un modelo 
bidimensional (Xle etaL, 1989). A menos que se indique explfcitamente otra cosa, 
en lo sucesivo utilizaremos el modelo bidimensional del sensor. Por lo tanto, las 
cargas electricas qt y las capacitancias cu y c,y de las ecuaciones (1) y (2) seran 
ahora cantldades por unidad de longltud de los electrodos en sentido axial. 
Utilizaremos la tilde (o sea, y c^) para denotar las cantidades totales que 

resultan al considerar la longitud real de los electrodos en la practlca. Las 
variables anteriores se relacionan entre sf por medio de la longltud de los 
electrodos L, de acuerdo con 



Si dividimos el interior del sensor bidimensional en p regiones de igual area 
(que llamaremos pixeles) donde la permitividad se considera constante, entonces 
la versidn discreta del problema directo es 
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^1 


c = 




= f (B) = 




con 8 = 





donde c es el vector de capacitancias mutuas (por unidad de longitud), ft son 
funciones no lineales que no se conocen explfcitamente y e es el vector de 
5 permitividades correspondiente a las p regiones o pixeles dentro de la zona de 
exploracion. 

Aplicando la Ley de Gauss, las capacitancias mutuas por unidad de longitud 
axial de los electrodos se pueden calcular como 

10 = - = ~^^(^V^')-^ = -^f^^^' (6) 

^ vj 7 2 V I dn 

donde So es una constante ffsica llamada permitividad del vacfo, igual a 
8.854 X 10"^^ faradios por metro, r, es una curva cerrada que rodea al electrode /, 
dl es un vector (normal) que representa un elemento infinitesimal de la curva Ft , 

15 dl es un elemento de longitud de diclia curva, el sfmbolo ' • ' representa el producto 
escalar de dos vectores, y ^ es el potencial electrostatico que se produce en el 
sensor al apllcar un voltaje de V voltlos en el electrode j (al que llamaremos 
electrode fuente o emisor) y 0 voltios en todos los demds (a los que llamaremos 
electrodos receptores). 

20 El potencial ^' esta determinado por la solucl6n de la siguiente ecuacidn 

diferencial parcial 

V*six,y)V^' =0 (7) 

25 sujeta a las condiciones de frontera (a) K voltios en el electrode fuente y (b) 
0 en los electrodos receptores y en la pantalla externa. En general, la 
ecuacion (7) no tiene solucion analitica y es necesario resolveria numericamente. 

Al problema de estimar cual es la distribucion espaclal de permitividad dentro 
del sensor que corresponde a un conjunto determinado de valores de 
30 capacitancia mutua, se le conoce como problema inverse, y es el problema que 
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deben resolver los m^todos de reconstmccion de imagenes. Normalmente, la 
estimacion de la permltividad se hace en forma discreta, representendola por un 
vector £ como el de la ecuacidn (5), que se debe calcular a partir de un vector de 
capacltancias mutuas observadas c, obtenido por medio de un aparato de 
medicion apropiado. 

Para resolver el problema in verso, la mayoria de los sistemas de TCE 
emplean el metodo de proyeccion Inversa lineal (LBP) (Plaskowski et al., 1995; 
Yang y Peng, 2003; Xie etaL, 1989, 1992), que se describe a continuacion. Como 
primer paso, para cada uno de los m = in(w-l) pares de electrodes que es poslble 
formar, se calcula un mapa de sensibilidad definido por 

s^^k) = ''^^^^^Jnen.p, pg^3 1, . . , n (8) 

donde k es el numero de pixel (de 1 a p), Ci(k) es la capacitancia medida con el par 
de electrodes / cuando el area correspond iente al pixel k esta llena de un material 
de alta permltividad y el resto del sensor contiene material de baja permltividad, 
mientras que cufuu) y cuemp) son las capacltancias del par de electrodes / cuando el 
Interior del sensor esta lleno de material de alta y baja permltividad, 
respectlvamente. Generaimente estos mapas de sensibilidad se calculan 
resolvlendo num§ricamente la ecuaci6n (7) y apllcando la ecuacion (6). 

Ya tenlendo los mapas de sensibilidad, 6stos se pueden usar para obtener 
una imagen de pennitividad a partir de cualquier vector de m = 
capacltancias mutuas medidas, c, conrespondientes a alguna distribuci6n 
desconocida de material dentro del sensor. Para ello se deben primero normalizar 
las lecturas de capacitancia medidas de acuerdo con 

;t^(;t) = (9) 

donde ^ es la capacitancia mutua normalizada conrespondiente al par de 
electrodes i y c/ es el valor de capacitancia medido en la practica para ese mismo 
par de electrodes. 
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La fdrmula bdsica del mStodo LBP calcula un 'nivel de gris' g(k) para cada 
pixel como 

III 

Sik) = para k^l,,. ,p (10) 

En principlo, se supone que ese nivel de gris esta relacionado linealmente 
con la pennitividad, con g I y g = 0 correspondiendo a la permitividad de los 
materiales de alta y baja penmitividad, respectivamente. El metodo LBP se basa 
en realizar una aproximacion lineal de un problema que, como ya se ha 
mencionado, en realidad es esencialmente no lineal (Gamlo y Ortlz-Aleman, 
2003), For lo tanto, este metodo de reconstruccion de imagenes provoca enrores 
considerables, los cuales son particuiarmente graves si las diferencias de 
permitividad son grandes. 

Hasta ahora, la principal alternativa al metodo LBP ha sido el empleo de 
metodos iteratlvos que buscan minimizar una funcion objetivo, usando tecnicas de 
optimizacion local como el metodo de Newton-Raphson regularizado o variantes 
del mismo (Yang y Peng, 2003). Como un ejempio de estos metodos esta el que 
se usa en el programa de compute EIT2D (Vauhkonen et al., 2001), desarrollado 
por investigadores de FInlandIa y el Reino Unldo. Su metodo de basa en 
minimizar, con respecto a e, el funcional 

llw-ccafcf + a'llLsIl' (11) 

donde a es un parametro de regularizacion, L es una matriz que contiene algun 
tipo de infonmacion a priori sobre la suavldad de c, y Ccaic= f(6) es el vector de las 
capacitancias mutuas calculadas para una distribucion de permitividad especifica 
dentro del sensor. En el programa EIT2D, Ccaic se calcula resolviendo el problema 
directo por medio del metodo de elementos finitbs. Partiendo de una estimaclon 
Iniclal 8o, la minimizaci6n se realize por medio del siguiente procedimiento iterative 
(bSsicamente un metodo de tipo Newton con regularizacion de Tikhonov) 

s*.t = 8, + [j/j, + a'l7l.y' {j; [c_- f(8,)] - a^L^Ls, } (12) 
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donde Jjt se conoce como la matriz jacobiana, que contlene las derivadas 
parciaies de f(8), evaluadas en Sk 



Sin embargo, estos metodos de reconstaicclon de Imdgenes tienen el 
problema de que requieren, para su buena operacion, de uno o m^s parametros 
de regularizacl6n cuyo valor correcto tiene una fuerte dependencia. precisamente, 
de las caracteristicas de la distribucion de permitividad que se pretende 
reconstruir, lo cual implicaria conocer de antemano la solucion del problema a 
resolver. Asimlsmo, generalmente estos metodos producen imagenes 
distorsionadas, debido a que la regularizacion tiene un efecto de suavizamiento 
excesivo sobre la permitividad obtenlda. Si la regularizacion es muy fuerte ocurrira 
el suavizamiento mencionado, y si es muy debil el metodo puede volverse 
inestables y/o no converger a la solucion buscada. 

Estos algoritmos de optimizacion local, durante su bOsqueda, exploran solo 
un pequeno sector del dominio de soluciones, restringido a la vecindad que 
circunda la solucion Inicial. SI la solucion optima del problema, es decir, el mfnimo 
absolute de la funcion objetlvo, se encuentra alejada de la solucion inicial, 
diffcilmente sera alcanzada debido a la presencia de minimos relatives 
interpuestos en su camino, lugares donde suelen quedar entrampados estos 
metodos, Los metodos mas utilizados correspondientes a esta categorfa son la 
inversion lineal por mmlmos cuadrados y las tecnicas que emplean el gradiente de 
la funcion objetivo, como el de maxima pendlente (en ingles 'steepest descent') y 
el gradiente conjugado. En general, los metodos de busqueda local explotan la 
escasa informacion derivada de la comparacion de una pequena cantidad de 
modelos, evitando asf una busqueda extenslva en el espacio de modelos 
(Sambridge y Drijkoningen, 1992). 

Los metodos de bOsqueda global, como su nombre lo indica, exploran todo 
el dominio de soluciones a lo largo del proceso de inversi6n. Hacen un rastreo 
exiiaustivo en el espacio de modelos. De esta manera, a pesar de existir 
soluciones parciaies del problema, es mayor la probabilidad de que la solucion 



J, = J(e,) = 




( i= l...m, 7=1.../? ) 



(13) 
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final corresponda al mejor ajuste entre los datos observados y los datos sinteticos. 
Esta clase de m^todos, en contraste con las t^cnicas locales, no requiere de la 
infomriaci6n proporclonada por las derlvadas de la funci6n objetivo debido a que 
en ellos el problema no requiere ser linealizado. Los algoritmos de optimlzacion 
global emplean criterios estocasticos para explorar simultdneamente todo el 
espacio de soluciones en la busqueda del modelo optimo. El m^s conocido de los 
metodos globales es el de Monte Carlo, que realiza una busqueda puramente 
aleatoria y no sesgada. Es declr, al generar cada nuevo modelo, no aprovecha la 
informacion obtenida de los modelos previamente evaluados (Gallagher et al., 
1991). La aleatoriedad no encauzada es el rasgo mas caracteristico de este 
metodo, que lo distingue del resto de los metodos globales. Entre las tecnicas de 
optlmizacl6n global, se encuentran tambien los metodos de algoritmos geneticos y 
de recocldo simulado. Ambos fueron concebidos como analogfas con sistemas de 
optimlzacion que ocurren en la naturaleza. Los algoritmos geneticos emulan los 
mecanismos de la evolucion bioiogica mientras que la tecnica de recocldo 
simulado tiene su base en la termodinamica. Ambos metodos son intrinsecamente 
no lineales y, por lo tanto, se prestan naturalmente a su aplicacion a la tomograffa 
de capacitancia, un problema no lineal. 

ESPECIFICACION DE LA INVENCION 

La presente invencion se refiere a un metodo de reconstruccion de imagenes 
para la visuallzacion de flujos multifasicos usando tomograffa de capacitancia 
el^ctrica, el cual esta basado en tecnicas heuristicas no lineales de optimlzacion 
global, en particular al metodo de recocldo simulado y los algoritmos geneticos. 
Aslmismo, considera el arreglo circular de electrodes metalicos rectangulares 
contiguos alrededor de la pared externa de un tube hecho de un material 
el6ctricamente aislante, formando un sensor A trav6s este sensor fluye una 
mezcia de fluidos en forma de flujo multifasico, cuya distribucion espacial en el 
interior del sensor se desea conocer. Se obtienen datos a base de efectuar 
mediciones de capacitancia mutua entre todos los electrodes de pared que es 
posible formar. Dichos datos dependen de la distribucion de los fluidos dentro de 
la tubena o tanque. 
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Por lo tanto, un objeto de la presente Invendon es proporcionar un metodo 
para procesar datos a base de efectuar mediciones de capacitancia mutua entre 
todos los electrodes de pared para reconstrulr una imagen de la distribucidn 
espacial de las fases o componentes de la mezcia multifasica, que fluye a traves 
del sensor, utillzando los m6todos de recocldo simulado y/o de algoritmos 
geneticos. 

Otro objeto de la presente invencion es que se utiliza en medldores 
tomograficos de flujo multifasico para cuantificar los diferentes fluidos (como gas, 
aceite y agua) producidos en un pozo petrolifero. 

Aun otro objeto mas de la presente invencion es que se usa para la 
visualizaclbn y monitoreo de otros flujos y procesos multifasicos que ocurren en la 
Industria, tales como, optimlzacibn del diseno y operacion de separadores, 
optlmizacidn del diseno de unldades de deslntegraci6n de fluldlzacion catalitica 
(FCC) y sistemas de lecho fluidlzado, y optlmizaci6n del diseno de reactores de 
lecho fijo. 
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BREVE DESCRiPCION DE LOS DIBUJOS DE LA INVENCION 

La figura 1 representa esquematicamente un sistema de tomograffa de 
procesos utilizado en la presente Invencion. 

La figura 2 representa esquematicamente un sistema de tomografia de 
capacitancia electrica empleado en la presente invencion. 

La figura 3 muestra un dibujo del sensor empleado en tomograffa de 
capacitancia usado en la presente invencion. 

La figura 4 muestra un corte transversal del sensor a la altura de la zona de 
medicion utilizado en la presente invencion. 

La figura 5 presenta un diagrama esquemdtico del metodo de recocido 
simulado empleado en la presente invencion. 

La figura 6 muestra el diagrama de flujo de la reallzacion del metodo de 
recocido simulado usado en la presente Invencion. 

La figura 7 presenta un diagrama esquematico del metodo de algoritmos 
genetlcos utilizado en la presente invencion. 

La figura 8 muestra esquemdtticamente el proceso de cruza de modelos 
usado en el mdtodo de algoritmos genetlcos empleado en la presente Invencion. 

La figura 9 muestra esquematicamente el proceso de mutaci6n de modelos 
usado en el metodo de algoritmos genetlcos en la presente invencion. 

La figura 10 muestra el diagrama de flujo de la realizacion del metodo de 
algoritmos geneticos utilizado en la presente invencion. 

La figura 11 muestra los resultados obtenidos al reconstruir imagenes de un 
flujo anular simulado de gas y aceite, usando el metodo de recocido simulado de 
la presente invencion. 

La figura 12 muestra los resultados obtenidos al reconstruir imagenes de un 
flujo estratificado simulado de gas y aceite, usando ei metodo de recocido 
simulado de la presente Invenci5n. 

La figura 13 muestra los resultados obtenidos ai reconstruir imagenes de un 
flujo simulado de aceite con burbujas de gas, usando el metodo de recocido 
simulado de la presente invencidn. 
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DESCRIPCION DETALLADA DE LA INVENCION 

Los nuevos procedimientos de reconstmccidn de imageries de la presente 
invencion estan basados en metodos heurfsticos de optimizacidn global, 
especfficamente el m§todo de recocido simulado (llamado en ingISs 'simulated 
annealing*) y los algoritmos gen6ticos. 

Una pluralidad de n electrodes metalicos se coloca alrededor de la periferia 
de una region a visualizar y se efectuan mediciones eiectricas de resistencia o de 
capacitancia entre ellos. O sea, se obtienen m = in(n-l) mediciones o dates. Las 
mediciones eiectricas son preferiblemente valores de capacitancia, aunque 
tambien podrian ser valores de resistencia. Como se muestra en la figura 3, los 
electrodes (B) pueden ser de forma rectangular y estar colocados sobre la pared 
externa de un tube (A) hecho de material elSctricamente aislante, formando asi un 
sensor, el cual contlene un fiujo multifasico o multicomponente (E). El objetivo es 
utilizar dichas mediciones para reconstruir una imagen de la distribucion espaclal 
del valor de la constante dielectrica (o permitividad) en la region a visualizar. 
Dicha distribucion de la permitividad refleja la distribucidn de las sustancias que 
ocupan el interior del sensor, donde se encuentra la regi6n a visualizar. 
Consideraremos que dicha region esta dividida en p partes iguales. 

Metodo de Recocido Simulado 

El metodo de recocido simulado esta basado en una analog la con el proceso 
termodinamico de la cristalizacion. Un fluldo mineral que se enfria lentamente 
hasta alcanzar un estado de energia bajo, da lugar a la formacion de cristales 
bien definidos. Si, por el contrario, la sustancia abandona su estado de equilibrio 
termico con un enfriamiento subito o parcial, entonces el cristal resultante tendra 
muchos defectos, en case de que la sustancia no forme un vidrio, caracterizado 
por el desorden metaestable de sus moleculas. Este concepto es utilizado en el 
contexto de los metodos de optimlzaci6n para reconocer configuraciones o 
modelos potencialmente dtlles. 

Los atomos de cada configuracidn molecular equivalen a los pardmetros del 
modelo en el problema inverse (o sea, la permitividad en los distintos pixeles de 
una imagen). La energfa del sistema para dicha configuracion se relaciona con la 
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funcion de costo (o de desajuste) asociada al conjunto de parametros que 
intervienen en el modelo. En nuestro caso para la presente invencion, la energia 
del sistema se asocia con la siguiente nomia Lz: 

m 

E^L^=^^ (i = l...,m) (14) 

donde f^meas son las m capacitancias mutuas medidas y c(Oc<i/c son las calculadas 
al resolver el problema directo para una distribucion de permitividad determinada 
8. Partlendo de una distribucion inicial de permitividad, el metodo genera una 
gama de configuraciones o combinaciones de parametros considerando una 
cierta temperatura T para el proceso. Para este proposito se emplea el criterio de 
Metropolis et aL (1953), que consiste en desplazar un parametro, en cada 
iteracion, una distancia aleatoria y pequena. Este desplazamiento provoca un 
cambio AE en la energia total del sistema. Si AE es menor o igual a cero, el 
desplazamiento del parametro es aceptado y la configuraclon resultante se 
considera como la nueva configuraclon actualizada. Cuando existe un incremento 
en la energia del sistema (AE es mayor que cero), la probabllldad de aceptacion o 
rechazo para el desplazamiento se determina como: 

P(AE) = 6-^/^ (15) 

Para saber si es o no admitido el cambio de posicion que implica el aumento 
de la energia del sistema, se elige en forma aleatoria un numero entre cero y uno, 
que se compara con el valor de la probabilidad correspondiente a AE. Si es menor 
dicho nOmero, se admlte el desplazamiento y se considera a la nueva 
configuraci6n como la actualizada; si es mayor, no se admlte y se conserva la 
configuraclon que se tenfa antes del movlmiento. Repitiendo sucesivamente este 
procedimiento se esta simulando el movlmiento t^rmico de los ^tomos de un 
sistema (que se encuentra en equilibrlo termico), a una temperatura fija 7. Si lo 
que se desea es alcanzar el estado base del sistema, es decir, el estado de 
menor energia y mayor ordenamlento, entonces se debe disminuir muy 
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lentamente la temperatura para simular un proceso cuasiestatico. Esto significa 
que, durante el enfriamiento, el sistema debe experimentar una sucesion de 
estados infinitesimalmente alejados del estado de equilibrio temnico. 

El m§todo de recocido simulado tiene tres componentes basicos 
(Vasudevan, 1991): una funci6n de energia o funcion de costo, una funcion de 
orden (criterio de Metropolis) y un par^metro que controla la temperatura del 
sistema. El proceso consta de tres ciclos anidados. 

En la figura 5 se muestra un diagrama de funcionamiento para el metodo en 
la presente invencion. El cicio externo (3) regula la temperatura del sistema. Cada 
vez que se cumple un ciclo, la temperatura del proceso disminuye al ser 
multiplicada por un factor Rr que normalmente es cercano a uno {0<Rt<1). De 
esta manera se lleva a cabo el enfriamiento lento y gradual que se propone. El 
cicIo intenmedio (2) se encarga de actualizar los valores, independientes entre sf, 
de una serie de constantes Ki asociadas a cada parametro. Dichas constantes 
determinan el mSximo incremento que podrS tener cada parametro al memento de 
ser perturbado en el cIcIo mas interne (1) del proceso. Los valores que adoptan 
estas constantes dependen de la cantidad de veces que liaya side aceptado el 
modelo actual al termino de cada secuencia de ciclos internes (1) segun el criterio 
de Metropolis. En el cicIo intemo (1) se perturban los valores de los parametros 
empleando los factores Ki , definidos en el ciclo intermedio (2). La perturbacion se 
realiza multiplicando a cada parametro por el resultado del producto de su 
correspondiente Ki por un valor aleatorio entre menos uno y uno. Posteriormente 
se calcula la respuesta sintetica del modelo actual y se evalua ei cambio de 
energia en el sistema asociado a la nueva configuracion de parametros. Dicha 
variacion de energia corresponde al desajuste entre la curve de dates sintetica y 
la observada o medida. Si ei desajuste decrece, entonces la nueva configuraci6n 
ser^ aceptada como la actual y perturbada de la misma manera. Si, por ei 
contrario, la perturbaci6n aleatoria produce un crecimiento en ei desajuste, 
asociado a un incremento en la energia E del sistema, a dicha configuraci6n se le 
asigna una probabiiidad de aceptacion de acuerdo con ei criterio de Metropolis. 

Los ciclos (1), (2), y (3) de la figura 5, se repiten, mientras la temperatura del 
proceso disminuye progresivamente. Conforme disminuye la temperatura, las 
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variaciones de los parametros son cada vez mas pequenas. De esta forma, la 
bQsqueda en el dominio de soluciones tiende a confinarse hacia los modelos 
asoclados con el minimo absolute de la funcldn de desajuste E. El resultado final 
es un conjunto de valores para los parametros (o sea, la permitividad en los 
distintos plxeles que forman la imagen) cuya respuesta sintetica reproduce los 
datos observados, salvo por un error suficlentemente pequeno. 

Como un ejempio, una realizacion especiTica, pero no la unica, del metodo 
de recocido simulado se presenta en la figura 6. En el bloque (1), partimos de una 
serie de valores iniclales de la temperatura, las constantes de perturbacion, las 
permitividades y la funcion de costo {To, Ki^o) V (con i = 1, y Eo). Al 

mismo tiempo, tambien en el bloque (1) se inicia el contador del cicio extemo, 
denotado por k. Posterionmente, en el bloque (2) se inician los contadores de los 
ciclos Interne e intermedio, i y j, respectlvamente, y entramos al cIcIo intemo. 
Como parte de dicho cicIo, en el bloque (3) se perturba aleatorlamente cada uno 
de los parametros (o permitividades), uno tras otro. Tambi6n en el bloque (3), 
cada vez que se perturba un parametro, se resuelve el problema directo y se 
calcula la funcl6n de costo o desajuste, E, aplicando la ecuaclon (14). Si hubo un 
decremento de E con respecto a la evaluacl6n anterior, se acepta el valor 
perturbado del parametro como el nuevo valor del mismo, se incrementa el 
contador de parametros / del cicIo intemo, y se precede a perturbar el siguiente 
parametro (si lo hay). Si, por el contrario, hubo un incremento de E, entonces se 
aplica el criterio de Metropolis, en el bloque (4), para decidir si se acepta o no el 
valor perturbado del parametro como su nuevo valor actualizado. Si, de acuerdo 
con dicho criterio, se acepta el nuevo valor, entonces se incrementa el contador 
del cicIo intemo i y se procede a perturbar el siguiente parametro (si lo hay). Si, de 
acuerdo con el criterio de Metropolis, no se acepta el valor perturbado, entonces 
no se Incrementa el contador i y se procede a perturbar nuevamente el mismo 
parametro. Una vez que se han actualizado todos los parametros Si, en el bloque 
(5) se ajusta el valor de las constantes Ki (que detemninan la forma en la cual se 
perturba a los parametros en el cicIo interne) y se incrementa el contador del cIcIo 
intermedio,/. JVk: determine el numero de veces que se repetira el cIcIo Intemo sin 
que haya una dlsminuclon de la temperatura. O sea, el ciclo intermedio consiste 
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en la repeticion del cicio Intemo Nk veces, pero con diferentes valores de Ki . Esto 
se hace con el fin de evitar que la temperatura descienda demaslado 
rapidamente, lo cual puede ser contraproducente en algunas apllcaclones del 
recocido simulado. Sin embargo, en la presente invencldn en particular no ocurre 
asi, y el valor de Nk se toma como uno. Al terminar el cicIo Intermedio, en el 
bloque (6) se reduce la temperatura como se indico en parrafos anteriores y se 
incrementa el contador del cicio extemo k. El procedimiento anterior complete se 
repite hasta que k alcanza el limite de iteraciones Nt, o antes si se alcanza un 
valor suficientemente bajo de la funcion de costo E. 

La presente invencion tambien se refiere a un metodo de reconstrucclon de 
imagenes basado en algoritmos geneticos, el cual se describe a continuacion. 

MStodo de Algoritmos Geneticos 

Los algoritmos geneticos, originalmente propuestos por Holland (1975), 
representan una evoluci6n del metodo de Monte Carlo para problemas de 
optlmizacion fuertemente no lineales. La busqueda del modelo 6ptimo se lleva a 
cabo explorando simultSneamente la totalidad del espacio de soluciones, 
empleando una regia de transicion probabilistica para guiar dicha busqueda. El 
proceso se Inicia a partir de un conjunto de modelos seleccionados 
aleatoriamente. 

Los parametros de cada modeio se transforman en codigo binario para 
formar cadenas denominadas cromosomas, sobre las cuales se aplican criterios 
de seleccion natural y genetica. Los procesos de seleccion, cruza y mutacion 
actualizan la poblacion de modelos, dando lugar a una nueva generacion de 
cromosomas, emulando la forma en que los sistemas blol6glcos evolucionan para 
producir organismos mejor adaptados al entomo. El proceso complete se repite 
hasta que la media de la funcion de ajuste se acerca al maximo ajuste para toda 
la poblaci6n. 

El diagrama de flujo de la figura 7 resume el proceso utiiizado para apllcar un 
esquema de Inversion basado en algoritmos gen§ticos, similar al descrito por 
Rodriguez-Zuniga et aL (1996) y Ortlz-Aleman etaL (2002, 2003) para la inversion 
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de parametros eldsticos del subsuelo a partir de datos de inclinacion del terreno y 
para la Inversl6n de datos aeromagnSticos, respectivamente. Los pasos basicos 
del proceso se descrlben brevemente a continuaclon, con referenda a la figura 7. 

Discretizacion 

En el bloque (1) de la figura 7, Ids parametros se representan por medio de 
un vector de incognitas s, que representa la distribucion de permitividad 
desconocida. La funclon de costo, que determina el ajuste entre los datos 
obsen/ados y la respuesta sintetica de un modelo dado, se denota como ^(s). La 
codlficacion de los parametros se realiza tomando en cuenta la extensi6n 
necesaria de la bOsqueda en el espacio de modelos y la resolucion deseada 
(Stoffa y Sen, 1991). De esta manera, la extensidn se define para cada parametro 
estableclendo un par de cotas at y fc/, es declr, ai<Si<bi, La resolucion se 
controla con el intervalo de discretizacion di , deflnldo como 



donde Ni es la cantidad de posibles valores para el parametro durante el proceso 
(Sambridge y Drijkoningen, 1992). Los modelos permitidos, 8, definidos por el 
conjunto de parametros €i , est^n restringidos al domlnio de valores 



Poblacion Inlcial 

Tambl§n en el bloque (1) de la figura 7, a partir de combinaciones aleatorias 
de los parametros, se construye una poblacion inicial de modelos, cuya dimension 
depende del problema particular a resolver. Cada combinacion se traduce en un 
conjunto de indices enteros, definidos por la ecuacion (17). Estos valores enteros, 
que establecen el valor particular de cada parametro del modelo, posteriormente 
se codifican como cadenas binarias llamadas cromosomas. Estas cadenas estan 
formadas por grupos consecutivos de bits, llamados genes, que representan el 
valor de los diferentes parametros Si . Los modelos que forman las poblaciones en 



(16) 



= a. + Jd^ para J = 0,..,,iV^ 



(17) 
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generaciones posteriores se crean a partir de los tres mecanismos evolutivos 
esenciales: la seleccion, la cmza y la mutacion. 

Problema directo y evaluacidn de la funclon de costo 

El problema directo, en el bloque (2) de la figura 7, conslste en calcular la 
respuesta teorica o sintetica para cada modelo en el proceso iteratlvo. 
Posteriomnente, en el bloque (3) de la figura 7, dicha respuesta sintetica se 
compara con las observaciones o dates mediante alguna medida de semejanza 
conocida como funclon de costo. El criterio utilizado en este trabajo es la norma La 
definida anteriormente en la ecuacion (14), aunque puede haber otros. La norma 
seleccionada para evaluar el ajuste debe ser sensible a la forma y complejidad de 
las curvas observada y calculada. 

Seleccion 

En el bloque (4) de la figura 7, a partir de una poblacion de Q individuos y de 
sus respectlvas funclones de costo E{Zk^ se determina para cada 

modelo una probabilidad ('acumulada') de seleccion P(zk) que dependera de su 
nivel de ajuste. 

Una formula que puede ser usada para determinar la probabilidad 
acumulada de seleccion es: 

donde E^ax, y Eprom son las funclones de costo maxima y promedio de la 
generacion, respectivamente, y g es el numero de individuos de la poblacion. 
Como siguiente paso, se utillza un procedimiento de ruleta sesgada (Goldberg, 
1989) para seleccionar una nueva poblacion de modelos. Se generan Q nCimeros 
aleatorlos ru entre cero y uno. Si P(sjt-i) P(8jt), entonces sjt es selecclonado 
para former parte de la nueva poblacion. Esto implica que en la nueva poblaci6n 
puede haber modelos 'gemelos'. Adicionalmente, se pueden ahadir a la poblacion 
de Q modelos, 'clones' de los mejores modelos, los cuales no se someteran a la 



23 



wo 2005/019779 



PCT/MX2003/000067 



cmza ni a la mutacion, con el objeto de asegurar que dichos modeios deseables 
no se pierdan al pasar por esos procesos que son de naturaleza aleatorla. 

Altematlvamente, la probabilidad de selecclon puede determlnarse como 
(Sambrldge y Drijkoningen, 1992) 

P(s,) = a^bEis,) (19) 
que describe una distribucion de probabilidades lineal, y 

P(s,) = Ae"^^^'"'^ (20) 

que corresponde a una distribucion exponencial. Los vaiores que suelen tomar las 
constantes a,b,AyB son los siguientes 

b = Q^' iE^ - E^^J-' , a > bE^ 

donde Emax, Eprom V Ecr son las funclones de costo maxima, promedio, y la 
desviacion estandar de todos los desajustes de la poblaclon Iniclal, 
respectivamente. 

Stoffa y Sen (1991) proponen un criterio de selecclon basado en una 
probabilidad de actuallzaclon. El criterio consiste en comparar el desajuste de 
cada modelo perteneciente a la generaci6n actual con el de un modelo de la 
generacion anterior, elegido al azar. Si el desajuste del modelo actual es menor, 
entonces §ste se conserva. De lo contrario, se considera un valor P„ que 
establece la probabilidad de sustituclon del modelo anterior por el actual. Este 
procedlmiento controla la Influencia de los ajustes de generaciones prevlas sobre 
la poblaclon actual. El valor sugerido por Stoffa y Sen (1991) para la probabilidad 
de sustituclon P„ es del 90%. 

Cruza 

Los nuevos modeios se engendran en el bloque (5) de la figura 7 a partir de 
una generacion progenitora de Q modeios. En forma aleatoria se integran QI2 
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parejas de individuos. Cada pareja es potencialmente capaz de cmzarse. Para 
determinar cuales de ellas llevaran a cabo la cruza, se le asocia a cada una un 
nOmero aleatorio entre cero y uno. Si dicho nOmero es menor que la probabilidad 
de caiza Pc entonces la pareja correspondlente efectOa el apareamlento. De lo 
contrario, la pareja se preserve intacta en la siguiente generacion. 

El mecanismo de cruza se basa en elegir aleatoriamente la posici6n de un 
gen para ambas cadenas. Las cadenas se parten en ese punto para intercambiar 
informacion entre ellas, tal como se observa en la figura 8, donde se muestra en 
(A) la particion de los cronnosomas padres en cuatro gametes, y en (B) el 
intercambio y union de dos gametes formando dos cigotos hijos. El proposito de 
cruzar dos cadenas diferentes es explorer nuevas regiones del dominio de 
soluciones, donde pudiera ubicarse el mfnimo absolute. En el proceso normal de 
cruza, las parejas que se aparean tienen dos hijos, y la poblaci6n de modelos se 
mantiene autom^ticamente en Q individuos. Aitemativamente, se puede 
considerar que las parejas tienen s6lo un hijo, en cuyo case la poblacion se debe 
completer con modelos generados aleatoriamente hasta volver a tener Q 
individuos. 

Mutacion 

La mutacion, al igual que la reproducclon sexual (o cruza), propicia la 
dlversidad genetica en una poblacion. La mutacion hace posible que la busqueda 
prospere cuando se encuentra confinada en los alrededores de un mfnimo local. 
La mutacion se realiza en el bloque (6) de la figura 7, mediante el cambio de 
paridad de un bit seleccionado al azar dentro de la cadena binaria (cromosoma). 
El porcentaje de modelos a los cuales se aplica el proceso de mutacion, al igual 
que en la cruza, depende de un parametro denominado probabilidad de mutacion, 
simbolizado por Pm . Este mecanismo previene la convergencia prematura del 
m6todo, cuando la poblaci6n es excesivamente homogenea e incapaz de 
continuar el proceso evolutive. En la figura 9 se representa una cadena arbitraria y 
se ejemplifica la mutacion del valor de un bit aleatorio, en (A) se muestra el 
cromosoma original y en (B) el cromosoma mutado. 
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Una altemativa para definir la probabilidad de mutaci6n fue propuesta por 
Yamanaka e Ishida (1996). Consiste en detenminar el nival de homogeneidad de 
los Indivfduos en cada generacion mediante el calculo de un coeficiente de 
vaiiacl6n promedio para cada pardmetro, mediante la expresibn 



donde p es la cantidad de parametros, es el promedio del /-esimo parametro. y 

<Ti es la desviacion estandar. A continuacion se define Pm como funcion de es 
decir 



donde Pua es la probabilidad de mutacion inicial. 

Con la mutaci6n concluye la secuencia de operaclones que define a un 
algoritmo genetico. DIcha secuencia se repite hasta satisfacer aiguna tolerancia 
preestablecida. 

Como ejempio, una reallzacidn especffica, pero no la unica, del metodo de 
algoritmos genetlcos se presenta en la figura 10. En el bloque (1), se parte de una 
poblacion de Q modelos generados aleatoriamente, cada uno representando una 
distribucion de permitividad. Al mismo tiempo, tambien en el bloque (1) se inicia el 
contador de iteraciones i. A seguir, se Inicia el contador de modelos k y comienza 
un cicio por medio del cual, para cada uno de los Q modelos, en el bloque (2) se 
resuelve el probiema directo y se calculan la funcidn de desajuste £(8a) y la 
probabilidad acumulada de seleccion P(8jO. Posteriormente, en el bloque (3) se 
realiza la seleccion de modelos con menor E, de acuerdo con el valor de PCs*) y 
aplicando el procedlmiento de ruleta sesgada ya descrito con anterioridad. Luego, 
en el bloque (4) se pasa a la cruza de modelos. En dicho bloque (4), se forman 
aleatoriamente Q/1 parejas, que se cruzan selectlvamente de acuerdo con el valor 
de Pc , produciendo solo un hijo por pareja. En ese mismo bloque (4), las parejas 
no cruzadas se conservan (ambos integrantes) iguales y la poblacion se completa 




(21) 



>^,^ para 7>0.1 

0.1 para 0.02 <;?^< 0.1 

0.2 para y < 0.02 



(22) 
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a Q individuos con nuevos modelos aleatorios. Por ultimo, en el bloque (5) se 
calculan el coeficiente de variacion promedio y la probabilidad de mutaclon P„, , 
y se reallza la mutacl6n de modelos de acuerdo con la descripci6n dada en 
pSn^afos anteriores. Todo el procedlmlento anterior se repite continuamente, 
incrementando i en cada iteraci6n. El algoritmo se temnina cuando se aloanza un 
numero preestablecido de iteraclones N o cuando la funcion de desajuste es lo 
suficientemente pequena. 

Fomriulacion del problema directo. 

El problema directo consiste en el calculo de las capacitancias mutuas cy , i ^ 
j\ que resultan de la presencia de una distribucion de permltividad 8 en el Interior 
del sensor. Tanto el metodo de recocldo simulado como el de algoritmos 
geneticos requieren de la solucl6n repetida del problema directo. Por lo anterior, 
es importante contar con un metodo adecuado para resolver dicho problema, que 
logre un equilibrio razonable entre exactitud (o precision) y rapidez. En el contexto 
de esta invencion, el problema directo se puede resolver utilizando una rutina 
optimizada desarrollada por los autores basada en el metodo de volumenes 
finitos. que es descrita muy brevemente a continuacion. Esta rutina es mas 
eficiente que las reportadas hasta hoy en la literatura (Yang y Peng, 2003) pues 
es comparable en precision con implantaciones basadas en el metodo de 
elementos finitos con mallas de 9,000 elementos triangulares. La velocidad de 
ejecucion es superior a las de los metodos de elementos finitos y diferencias 
finitas. La rutina est^ escrita en Fortran 90 y es totalmente transportable (se ha 
probado en computadoras tlpo PC y "clusters" de PCs basados en Windows y 
Linux, en estaclones de trabajo tipo SUN y Alpha y en supercomputadoras Cray). 
La rutina puede extenderse al case tridimensional sin mayores modiflcaciones y 
es Mcilmente paralelizable. 

El problema directo se resuelve mediante el m§todo de volumenes finitos en 
una configuracidn cilfndrica. De esta manera se elimlnan las soluclones 
indeterminadas en el centre del disco (como sucede en el metodo de diferencias 
finitas) y el refinamiento de la malla se vuelve mas flexible en comparaci6n con los 



27 



wo 2005/019779 



PCT/MX2003/000067 



metodos de elementos finitos. Se resuelven la ecuaci6n (7), que se repite a 
continuacion 

V-fir(jc,3;)V^*=0 

donde s es la permitividad y es el potenclal electrostatico que se genera 
cuando el electrode k es el emisor (o fuente). La ecuacl6n esta sujeta a las 
condlciones de frontera (a) j**=F voltlos en el electrode fuente y (/>) = 0 en los 
electrodos receptores y en la pantalla externa. 

Deflniendo las coordenadas radial y angular come r y 0, y utilizando el 
metodo de volunnenes finitos, la ecuacion discreta es formulada en forma 
conservadora en cada celda Slg como 

V.(^V(^*) Ji^^ = 0 para /=l...,iVr y (23) 

donde / y j se refieren a la discretlzacion en r y ft respectivamente, y NrV Ne son 
el numero de secciones en las que se subdivide el radio y la circunferencia del 
sensor, respectivamente. 

Aplicando el teorema de Gauss en coordenadas polares, las ecuaciones 
discretas pueden escribirse como 

j^^eV^'^dr^j = 0 (24) 

donde r^- es la frontera de la celda de volumen finito Cly . La frontera Fij se define 
por r^F y a lo largo de las coordenadas radiales, y Tat y r^- a lo largo de las 
coordenadas angulares. La ecuacion (24) puede expresarse como la suma de los 
flujos a traves de las caras Tj^ , , y Ffr 
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De (25), el termino correspondiente a los flujos en el radio cero desaparece y 
el problema es equivalente a resolver las ecuaciones en la cercania del centre 
sobre triangulos que tienen un vertice en el centre. Entonces, el sistema discrete 
de ecuaciones para el problema directo resulta bien planteado. El sistema 
complete es similar a un sistema de ecuaciones Laplaciano y debe resolverse un 
sistema diagonal bandeado que incluye las condiciones de frontera peri6dicas 
impuestas per la geometrfa del problema. La matriz correspondiente es positiva 
definlda y no simetrica, caractensticas que se pueden aprovechar para 
selecclonar especfficamente los metodos optimos de solucl6n. 

Finalmente, las capacitancias mutuas se calculan integrando los gradientes 
del potencial a lo largo de una curva cerrada que rodea los electrodes, de acuerdo 
con la ecuacion (6), que se repite a continuacion 




donde So es la permitividad del vacio (8.854 x lO'^^ faradios por metro), Ti es una 
curva cerrada que rodea al electrodo z, dl es un vector (normal) que representa un 
element© infinitesimal de la curva T/, dl es un elemento de longitud de dicha curva 
y es el potencial electrost^tico que se produce en el sensor al aplicar un voitaje 
de V voltios en el electrodo k (fuente) y cero voltios en todos los demds 
(receptores). La integracion se realiza utilizando una regia trapezoidal y los 
gradientes de potencial se calculan al cuarto orden. 

Durante el procedimiento de reconstrucci6n de una imagen de permitividad 
por el metodo de recocido simulado, es necesario resolver ei problema directo y 
encontrar el potencial electrico repetidamente para sucesivas distribuciones de 
permitividad relativamente simllares, mientras el m§todo converge a la solucion 
final. En virtud de que el potencial correspondiente a dichas distribuciones 
sucesivas cambia relativamente poco, es posible acelerar todo el proceso si se 
emplea un metodo iterative para resolver el problema directo, tomando como 
estimacion inlcial del potencial a la solucion correspondiente a la configuracion de 
permitividad anterior. Como ia estimacion inicial del potencial estara muy cerca de 
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la solucion, dicho metodo iterativo convergera en pocas iteraciones, alcanzando 
rdipidamente una exactitud aceptable. 

EJEMPLO 

Con objeto de evaluar el desempeno de los metodos de reconstruccl6n de 
imagenes descritos anteriormente, calculamos un conjunto de datos TCE 
sinteticos con ia subrutina del problema directo. Para emular las fuentes 
principales de incertidumbre en los datos producidos por uri instrnmento de 
medicion operando en condiciones de trabajo (esto es, los errores aleatorlos 
inherentes al proceso de medicion de cualquier variable fisica y los errores 
debidos a la precision propia del instrumento de medicion), los calculos de las 
capacitancias sinteticas para el modelo ideal fueron evaluadas con una precision 
numerica del orden de 10'^'^ en el metodo iterativo utillzado en el problema directo 
para el calculo de los potenclales. Con objeto de simular la Impreclsi6n 
(sistemStica) asociada al sensor TCE, durante el proceso de inversion (estimacion 
de la distribucion de permltlvidad elSctrIca en el interior del tubo) se utilizo una 
precision sensiblemente menor (10"^) en el calculo de los potenclales en el 
problema directo. 

Las interpretaciones sin restricclones de datos de campos potenclales 
pueden ser de muy poco interes practice debido a la gran ambiguedad presente 
entre las observaciones y las soluciones estimadas. La no unicidad en los 
problemas con campos potenclales surge principalmente de dos fuentes: la 
primera es la ambiguedad inherente provocada por la fisica del problema que 
permite la existencia de una inflnldad de soluciones que reproducen la anomalfa 
de campo potencial; la segunda resulta de la utillzacidn de un numero finite de 
datos que estan contamlnados por errores y que pueden contener lnformaci6n 
insuficiente para construir una solucidn Cinica del problema. Las estrategias que 
permiten superar esta no unicidad consisten en la incorporacion de suficlente 
informaci6n a priori para restringlr las soluciones resultantes a una regl6n del 
espacio de parametros que es conslderada ffslcamente razonable (Pllkington, 
1997). En el case particular de la tomograffa de capacitancia electrica se cuenta 
con Informacion acerca de los valores tfplcos de la permitividad electrica para el 
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gas, aceite y agua. El conocimiento de estos valores, asf como el contraste 
signlficativo entre las propiedades del gas y el aceite, permiten apllcar las tecnicas 
de inversion que aquf se discuten a dates de TCE para flujos bifaslcos 
contaminados por errores relatives menores al 2% (siendo este el nivel maximo 
de error que producen los sistemas de adquisicion de dates empleados 
normalmente en TCE). Por otra parte, si se cuenta con una estimacion estadistica 
precise de las incertidumbres en los datos, entonces es posible construir un 
nnodelo que considere las incertidumbres en los datos y en las estimaciones de 
los parametros, utillzando el esquema que proponen diversos autores para la 
inversion de datos de campos potenciales (por ejempio, Sen y Stoffa, 1995). 

Usando la subrutina del problema directo a base de volumenes finitos, se 
calcularon las capacitancias sinteticas considerando tres distribuciones 
especificas de permitividad, correspondientes a tres patrones de flujo tfpicos 
(anular, estratificado y de burbujas). Simuiamos un sensor TCE de 12 electrodes y 
los valores de capacitancia para todas las comblnaciones posibles de electrodes 
fueron calculados. Esos fueron los datos simulados. Consideramos una 
distribucion de dos componentes con un material de permitividad baja de 1 (gas) y 
otro de permitividad alta de 2.5 (aceite). Nuestras pruebas numericas se han 
hecho con datos libres de ruido y con datos contaminados con errores aleatorios 
de hasta 1%. El algoritmo esta codificado en Fortran 90 y conre en una 
computadora Pentium 4 de 1.7 gigahertz con 512 megabytes de memoria RAM. 
Probamos con una malia de 120^60 elementos para reducir los tiempos de la 
inversion (^ 30 minutes para 60,000 iteraciones), pero los resultados son validos 
para dimensiones mas grandes. La validez de los resultados de estas 
slmulaciones fue corroborada con resultados experimentales de laboratorio 
obtenidos empleando modelos fislcos y un sistema real de tomografia de 
capacitancia. 

Despues de una parametrizacion adecuada, tanto el metodo de algoritmos 
gen6ticos como el de recocldo simulado produjeron resultados muy similares, 
para los tres patrones de flujo conslderados. La calidad de las imagenes de 
permitividad reconstruidas depende principalmente del numero de iteraciones del 
m6todo, tal y como sucede para muchas otras aplicaciones (por ejempio, Ortiz- 
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Alem^n ef a/.. 1999. 2001, 2002, 2003; Cruz-Atienza, 1999). En la figura 11 se 
presenta la imagen ideal o Verdadera' de un flujo anular simplificado (A), y las 
reconstnjcciones obtenidas despuSs de 30,000 (B) y 60,000 (C) calculos del 
problema directo (o iteraciones del metodo). En las figuras 12 y 13 se presentan 
figuras similares para un flujo estratificado y un flujo de burbujas, 
respectivamente, mostrando en (A) la Imagen Ideal, en (B) la Imagen reconstoilda 
con 30,000 Iteraciones y en (C) la Imagen reconstrulda con 60,000 Iteraciones. En 
todos los cases mostrados en las figuras se uso el metodo de recocido simulado 
(con los algoritmos geneticos no hubo diferencias significativas en los resultados) 
y se partio de una distribucion de permitividad inicial homogenea. Estos 
resultados muestran claramente que las imagenes reconstruldas son muy 
parecldas a las Imagenes ideales o Verdaderas'. La exactitud de estas 
reconstrucciones es conslderablemente mejor que las reportadas hasta ahora en 
la llteratura (Yang y Peng, 2003). Si blen los metodos de esta Invencl6n no 
requieren, para converger, que la estimacidn inicial de la distribuci6n de 
permitividad este pr6xlma a la solucldn, es posible acelerar un poco el proceso si 
se emplea como estimacion inicial la imagen resultante de un metodo directo 
simple como el LBP. 
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NOVEDAD DE LA INVENCI6n 

Habiendo descrito la presente invenci6n, esta se considera como novedad y 
por lo tanto, se reclama como propiedad el contenido en las siguientes cidusulas: 

1 Un metodo de reconstrucclon de imagenes que comprende: (a) la obtencl6n 
de dates medidos correspondientes a la imagen; (b) el procesamiento de 
dichos datos medidos por medio de metodos heuristicos no lineales de 
optimizacion global, especificamente el metodo de recocido simulado o los 
algoritmos geneticos, para obtener una Imagen; y (c) el desplegado de dicha 
imagen en un dispositive de visualizacion. 

2. - Un metodo de conformidad con la clausula 1 , caracterizado porque los datos 

medidos son parametros el^ctricos en el perimetro de una region tal como el 
interior de una tuberfa, el de un poze o el de un tanque, para obtener una 
imagen de la distribucion espacial de la penmitividad el6ctrica (o constante 
dlel§ctrica) o de la conductividad electrica en la que se refleja la distribucion 
espacial de los materiales o sustancias, tales como gases y/o ifquldos. 

3. - Un metodo de conformidad con la clausula 2, caracterizado porque los 

parametros electrlcos son valores de capacltancia electrica medidos entre 
los distintes electrodes un sensor formado por una pluralidad de electrodes 
celocades a lo largo del perfmetro de la region (tuberfa, poze, tanque). 

4. - Un metodo de conformidad con las clausulas 2 y 3, caracterizado porque la 

Imagen obtenida es una distribucion espacial de la permitividad electrica (e 
constante dielectrica) en la region de inferos (tuberfa, poze, tanque), la cual 
refleja la distribucidn espacial de los materiales o sustancias, tales como 
gases y/o Ifquldos, que ocupan dicha regidn (tuberfa, poze, tanque). 

5. - Un metodo de conformidad con las clausulas 2 a 4, caracterizado porque la 

imagen esta formada por una cantidad finita de subregiones o pixeles en la 
visualizacion, cuyo nOmero depende de la resoluci6n deseada. 

6. - Un metodo de conformidad con la clausula 3, caracterizado porque el sensor 

esta formado por un tubo fabricado de un material electricamente alslante, 
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sobre cuya pared externa se coloca un anillo de electrodes met^Hcos 
rectangulares. 

7. - Un metodo de conformidad con las cl^usulas 3 y 6, caracterizado porque el 

sensor contiene un flujo multifaslco o multicomponente, y obtiene Im^genes 
que muestran la dlstribuci6n de las fases o componentes, tales como gases 
y/o Ifquidos de acuerdo a lo dicho en antecedentes. 

8. - Un metodo de conformidad con la clausula 1 , caracterizado porque empiea 

el metodo de recocido simulado como el metodo de optimizacion. 

9. - Un metodo de conformidad con la clausula 8, caracterizado porque se 

minimiza iterativamente una funcion de costo asociada con la energia del 
sistema, con respecto a e (el vector de los valores de la permitividad en cada 
pixel de la imagen), y siendo de la forma: 

= = — — nrz 1, . . , m) 

donde c^^ son las m capacltanclas mutuas medldas y cf^Csj^) son las 
calculadas al resolver el problema directo para una distribucion de 
permitividad determinada e^. 

10. - Un metodo de confonnidad con la clausula 9, caracterizado porque para 

efectuar (a minimlzacion se empiea el criterio de Metropolis. 

11. - Un metodo de confomiidad con la clausula 10, caracterizado porque en el 

proceso iterative de minlmizacl6n se empiea como estimacion inicial de la 
distribuci6n de permitividad, 8, tanto una distrlbucl6n liomogdnea como la 
que resulta de aplicar a los dates medidos el metodo de proyeccion inversa 
lineal (LBP). 

12. - Un metodo de conformidad con la clausula 9, caracterizado porque el calculo 

de las capacitancias calculadas cf^^(Sk) , conocido como el problema directo, 

se realiza por medio del metodo de volumenes finitos. 

13. - Un metodo de conformidad con la clausula 12, caracterizado porque para la 

solucion del sistema de ecuaciones que resulta al resolver el problema 
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directo (cSlculo de cf^'^Ce^)), se utilizan m^todos iterativos que convergen 

rapidamente al emplear como estimaci6n inicial del potencial electrostdtico, 
el potencial que resuito de la solucion del problema directo en la iteracldn 
anterior del problema inverso. 

14. - Un metodo, de conformidad con la cl^iusula 1, caracterizado porque se 

emplea el metodo de algoritmos geneticos como el metodo de optimizacion. 

15. - Un metodo de conformidad con la clausula 14, caracterizado porque, a partir 

de una poblaclon inicial de Q modelos de permitividad e^o} 
apllcando mecanismos evolutivos como seleccion, cruza y mutacion, se 
obtienen Iterativamente nuevas poblaciones. 

16. - Un mdtodo de confomiidad con la clausula 15, caracterizado porque los 

indivlduos poseen una funcion de costo o de desajuste pequena, donde la 
funci6n de desajuste est& dada por: 

donde c"'^"" son las m capacitancias mutuas medidas y cf^\i&j^) son las 
calculadas al resolver el problema directo para el modelo (o distribucion de 
permitividad) Zk- 

17. - Un metodo de acuerdo con las clSusulas 15 y 16, caracterizado porque la 

probabilidad acumulada de seleccion para un modelo determinado Zk esta 
dada por 



donde Emax^ y Eprom son las funclones de costo maxima y promedio de la 
generacion, respectivamente, y g es el nOmero de indivlduos de la 
poblacion. 
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18. - Un metodo de conformldad con la clausula 17, caracterlzado porque se 

puede emplear el procedimiento de la ruleta sesgada para decidir cuales 
modelos son selecclonados en cada iteracion. 

19. - Un metodo de confomiidad con la clausula 18, caracterizado porque la caiza 

y la niutaci6n de modelos se realizan aleatoriamente de acuerdo las 
probabilidades de cruza y de mutacion, Pc y Pm - 

20. - Un metodo de conformidad con la clausula 19, caracterizado porque la 

probabilidad de mutacion Pm esta determinada mediante el coeficlente de 
variacion promedio dado por 



donde p es la cantidad de parametros, ^. es el promedio del /-esimo 

par^metro, y Oi es la desviaci6n estandar. 
21.- Un metodo de conformldad con las clausulas 19 y 20, caracterizado porque 
Pm se define como funclon de es decir 



donde Pm es la probabilidad de mutacion inicial. 
22.- Un metodo de conformidad con las clausulas 15 y 16, caracterizado porque 
el calculo de cf^'iZj^), conocldo como el problema directo, se realiza por 
medio del metodo de volOmenes finltos. 




>^,^ para y>OA 

0.1 para 0.02 < / < 0.1 

0.2 para y < 0.02 
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